Multifractal nature of the surface local density of states in three-dimensional 
topological insulators with magnetic and nonmagnetic disorder 
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We compute the multifractal spectra associated to local density of states (LDOS) fluctuations 
due to weak quenched disorder, for a single Dirac fermion in two spatial dimensions. Our results 
are relevant to the surfaces of Z2 topological insulators such as Bi2Se3 and Bi2Te3, where LDOS 
modulations can be directly probed via scanning tunneling microscopy. We find a qualitative differ- 
ence in spectra obtained for magnetic versus non-magnetic disorder. Randomly polarized magnetic 
impurities induce quadratic multifractality at first order in the impurity density; by contrast, no 
operator exhibits multifractal scaling at this order for a non-magnetic impurity profile. For the 
time-reversal invariant case, we compute the first non-trivial multifractal correction, which appears 
at two loops (impurity density squared). We discuss spectral enhancement approaching the Dirac 
point due to renormalization, and we survey known results for the opposite limit of strong disorder. 

PACS numbers: 73.20.-r, 73.20. Jc, 64.60. al, 72.15.Rn 



I. INTRODUCTION 

The defining attribute of a 3D Z2 topological insulator 1 
(TI) is the presence of an odd number of 2D massless 
Dirac bands at the material surfaced Unlike the Dirac 
electrons that can appear in a purely 2D system (no- 
tably in graphene), the surface states of a (strong) 3D TI 
are robustly protected from the opening of gap, so long 
as time-reversal symmetry is preserved. The protection 
can be viewed as a consequence of the parity anomaly^— 
which "holographically" links surface states separated by 
a topologically non-trivial bulk, and gives rise to the sig- 
nature properties of the Z2 TI state: the half-integer 
quantum Hall effect, quantized magnetoelectric coupling, 
"axion" electrodynamics, etci2a2 As stressed by Schny- 
der et al. in Ref. 0, the robust character of the surface 
states in the presence of quenched disorder can also be 
taken as a principal characteristic of a topological insu- 
lator. In particular, these states are protected from An- 
derson localization,— even in the presence of a "strong" 
impurity potential, so long as time-reversal invariance is 
preserved i^iiS 

With its 2D Dirac band pinned to an exposed surface, 
a 3D TI is ideally suited to local probes such as scanning 
tunneling microscopy (STM). In spectroscopic mode, an 
STM captures an areal map of the local density of states 
(LDOS). There are several ways of analyzing such data. 
One is to look for quasiparticle interference (QPI)^i~— in 
the LDOS Fourier transform. This method is useful for 
determining short-distance details, and contains similar 
information as an analysis of LDOS Fricdcl oscillations 
in the presence of a single impurity^ It has been applied 
in TIs to experimental data and analyzed theoretically in 
Refs. [H[ll and [Hill, respectively. In QPI, the disorder 
is employed primarily as a facilitator to gleam informa- 
tion about the clean system^ 

Multifractal analysis^— provides a complementary 



method better suited to extracting large-distance, 
disorder-dominated features in the same LDOS data 
field. It is a standard tool for assaying quantum interfer- 
ence phenomena, and is employed in the analysis of wave- 
functions near a metal-insulator transition^— as well as 
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FIG. 1: Sketch of disorder "flavors" on the surface of a 
Z2 topological insulator. In the time-reversal invariant case, 
the impurities are neutral adatoms or charged dopant ions, 
depicted as spheres in (a) . The effects of these on the surface 
Dirac theory [Eq. (|3.1[) ] are encoded in the scalar potential 
V(r). In the case of magnetic disorder, the impurity spins 
are indicated by the arrows in (b) and (c). In the limit that 
the spins reside in the plane of the surface, (b), the disorder 
appears as a vector potential A(r). The opposite case of out- 
of-plane polarization, (c), gives the random mass M(r). The 
case of generic time-reversal breaking disorder has all three 
potentials present. 
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mesoscopic fluctuations in diffusive metallic systems i 22 ' 23 
In this paper, we derive new results for LDOS multifrac- 
tal spectra associated to disordered topological insulator 
surface states. In particular, we extend the pioneering re- 
sults of Ref. [24] to the generic cases of time-reversal (T) 
preserving and breaking impurities. Our calculations are 
performed in the near-ballistic limit^ wherein weak dis- 
order enters as a perturbation to the clean Dirac band 
structure. A key characteristic of 2D Dirac fcrmions is 
that this weak disorder regime is continuously connected 
to more conventional domains of multifractal analysis, 
i.e. the diffusive (symplectic) meta l 20 ' 23 and the inte- 
ger quantum Hall plateau transition! 18 ' 26 " — These ap- 
pear at strong coupling (many impurities) for dirty Dirac 
fermions ,£ i 10 i 24 i 29 

We consider the case of a single flavor Dirac surface 
band, relevant to (e.g.) the TIs Bi 2 Se3 and Bi 2 Te 3 
(Refs. iim. The different kinds of T -preserving and 
T-breaking disorder are sketched in Fig. Q] We demon- 
strate that the LDOS multifractal spectra observed in 
the absence of time-reversal symmetry breaking (i.e., for 
non-magnetic disorder) is qualitatively weaker than that 
induced by magnetic impurities. In particular, the first 
multifractal correction obtains at first order in the im- 
purity density for the case of broken T, while the first 
non-trivial amplitude appears at second order in the T- 
invariant case. We compute the leading terms via one- 
and two-loop calculations, respectively. We also compute 
unnormalized spectra for the spin LDOS 3 ^ in the case of 
magnetic impurities. We show that renormalization ef- 
fects can enhance multifractality near the Dirac point. 
Finally, we summarize prior results on various strong- 
coupling regimes. Our goal is to sketch the full portrait 
of quantum interference physics on the surface of a TI, 
valid when intcrparticle interactions can be neglected. 

Our results indicate that the long-distance, disorder- 
dominated features captured by the multifractal analysis 
behave in many cases opposite to the short-distance char- 
acteristics that appear in quasiparticle interference*^— 
In Ref. 0, the authors observed that QPI is strongest 
for the spin LDOS response to magnetic impurities, while 
the unpolarized LDOS pattern vanishes for magnetic dis- 
order (in the first Born approximation). The QPI re- 
sponse of the LDOS to non-magnetic disorder is weak 
but non-zero.— By contrast, in this work we find that 
the LDOS multifractality is strongest for magnetic im- 
purities, while the spin LDOS spectrum comparatively 
exhibits the same or weaker strength fluctuations, de- 
pending upon the polarization direction. 

The weak influence of non-magnetic disorder is tied to 
the intrinsic spin-orbit coupling that defines the mass- 
lcss Dirac kinetic term. Multifractality is suppressed at 
one loop due to interference mediated by the Dirac pseu- 
dospin, which is proportional to the physical spin on a 
Z 2 insulator surface. The spin is also responsible for the 
suppression of backscattering from a single non-magnetic 
impurity.— On the TI surface, magnetic disorder Zeeman 
couples directly to the Dirac spin, enabling backscatter- 



ing in near-ballistic transport, and inducing multifractal 
LDOS fluctuations at the lowest order in the impurity 
density. 

A notable problem in experiments probing topological 
insulator surface states has been the unintentional dop- 
ing of carriers into the bulk bands, which then dominate 
transport measurements in large samples.— Even if the 
chemical potential is moved into the gap, it may reside 
far from the Dirac point, making it difficult to observe 
surface state carrier dynamics at low densities. In this 
respect, STM offers several advantages over transport ex- 
periments. First, the position of the chemical potential is 
no barrier to probing states at the Dirac point, since the 
latter can always be reached by tuning the bias voltage 
(although the Dirac point is not guaranteed to reside in 
the bulk gap)* 3 ^ Assuming that the Dirac point or the 
low density regime can be accessed by tuning the tunnel- 
ing bias, the advent of a finite, even large doping of the 
surface and/or bulk states may actually play a beneficial 
role in facilitating the observation of disorder-induced 
quantum interference effects. This is because a finite 
carrier density screens the long-range Coulomb potential 
introduced by charged defects. The potential landscape 
formed by screened impurities is short-range correlated 
on scales larger than the screening length. Good screen- 
ing eliminates the problem of electron and hole puddle 
formation ) 34 ' 35 which has until recently 3 ^ occluded trans- 
port and other properties of Dirac carriers in graphenc 
near the Dirac point. On the other hand, a low density of 
poorly-screened bulk dopants induces a long-range corre- 
lated potential and puddle formation, as in grapheneJ^ 
LDOS fluctuations in the puddle regime are an important 
topic for future work. 

Three-dimensional topological insulators provide us 
with an interesting paradigm flip for quantum interfer- 
ence phenomena. Isolating the surface state contribution 
in transport measurements is problematic. By compari- 
son, direct LDOS imaging is easier than in conventional 
semiconductor systems, wherein the 2D electron gas is 
typically buried in a layered material stack. Moreover, 
the amount of surface disorder can to some extent be 
controlled; for example, magnetic impurities can be de- 
posited across the surface of an otherwise high-quality 
bulk 3D sample. These can be charge-neutral adatoms 
or charged dopants; an example of the former (latter) is 
provided by iron (manganese) 3 ^ in Bi2Se3. 

This paper is organized as follows. We begin in Sec. [TT] 
with a lightning review of multifractal composite and spin 
LDOS measures. In Sec. IIII| we present new results for 
multifractal LDOS fluctuations in TI surface states, in 
the presence of weak disorder. We also show how renor- 
malization can enhance multifractality close to the Dirac 
point. Finally, in Sec. IIV1 we review previous results on 
various strong disorder regimes relevant to the I2 TI sur- 
face states and LDOS statistics. In particular, we discuss 
the symplectic metal, the integer quantum Hall plateau 
transition, and the Anderson insulator. Various technical 
details are relegated to appendices. In Appendix \^ we 
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review the symmetry classes of Anderson (dc) localization 
that appear in the disordered Dirac surface theory. In 
Appendix [B] we supply some details of our perturbativc 
calculations. 



II. MULTIFRACTAL LDOS MEASURES 

A. Definitions 

We suppose that the tunneling local density of states 
(LDOS) v{e, r) is imaged at a fixed energy e over anlxi 
field of view. The field is finely partitioned into a grid of 
boxes. The box edge length a <C L must be chosen larger 
than any "microscopic" scale l m , such as the correlation 
length of the random potential^ One introduces the box 
probability 



/ d 2 rv(e,r) 



A»n(e) 



A n 



J d 2 rv(e,r) 

A, 



(2.1) 



where A n denotes the n th box. LDOS multifractality 

is defined through the inverse of the participation ratio 
(IPR)^ 



V q {e) 



^)~GD T(9 ' £) . (2-2) 



The right-hand side (scaling limit) obtains when l m <C 
a <C L; corrections are down by higher powers of 
a/L. The exponent r(g, e) is the multifractal moment 
spectru m 18 ' 38 for LDOS fluctuations at energy e. 

The construction in Eqs. (|2.ip and (|2.2[) is useful for 
characterizing a system with extended states, or for an 
Anderson localized system in which L -C 6oc(e); £ioc de- 
notes the localization length. In what follows, we as- 
sume experiments are performed at sufficiently low tem- 
peratures so that inelastic cutoffs to quantum interfer- 
ence can be ignored^ 3 ^ A clean system with plane wave 
states at energy e has r(q,e) = 2(q — 1). Multifractality 
refers to the incorporation of corrections non-linear in q. 
Physically, these arise due to quantum interference via 
multiple scattering of electron waves in a dirty environ- 
ment, processes that serve as the precursor to Anderson 
localization,!^^ 

For weak disorder, the spectrum is typically dominated 
by the quadratic correctio n 18 ! 20 



T(q,e) = 2(q-l)-6(e)q(q-l), 



(2.3) 



where 9 > gives a measure of the disorder strength. As 
an example, in a weakly disordered 2D metal (with L <C 
£i oc for the orthogonal or unitary classes), one find a 20 ' 21 ' 23 



r 1 



2ir 2 N(e)D' 



(2.4) 



where N(e) denotes the average density of states, D is 
the classical (Drude) diffusion constant, and /3 £ {1,2,4}, 
depending upon the presence or absence of time-reversal 
symmetry and spin-orbit scattering] 23 ! 40 At stronger dis- 
order, higher order corrections in 9q must be retained; for 
the diffusive metals, results are known to four loops.— 

An alternative characterization of LDOS multifractal- 
ity is provided by the singularity spectru m 18 ' 38 /(a): 
Over a subset of the sample grid area that scales as 
(L/a) fi - a \ the box probability \i ~ (a/L) a . The sin- 
gularity spectrum is the Legcndrc transform of r(q). 



/(a) = qa - r(q), 



dr(g) 
dq 



For the quadratic spectrum in Eq. (|2.3|) . one obtains 

(2.5) 



In this "parabolic approximation," the strength of 
the multifractality is encoded in the peak position ao 
[f( a o) = 2], and the width aw of the spectrum such that 
f{a ± aw/2) = 0, 



a =2 



a w = iVW. 



(2-6) 



Part of the power of multifractal analysis for disordered 
quantum systems derives from the fact that the spectra 
[t(<?) or /(a)] typically depend only upon a few gross 
measures of the impurity potential. In the case of dirty 
metals, the entire spectrum can be computed as an ex- 
pansion in one parameter, the inverse conductance (con- 
sistent with scaling theory) i 8 ' 20 ' 22 At a non- interacting 
Anderson localization transition, r(q) and /(a) become 
universal functions, so that the critical point is charac- 
terized by an infinite set of critical exponents [e.g., the 
expansion coefficients for r(q)]. 

The spectra above have been defined for data collected 
in a single fixed realization of the disorder. Strictly 
speaking, Eq. (|2.3p then applies only for \q\ < q c , where 
q c = \j2jQ. Outside of this range, the r(q) associated 
to a fixed disorder realization is linear, a phenomenon 
known as spectral termination^ - — [This assumes that 
higher order corrections can be ignored for q > q c . Re- 
gardless, the r(q) spectrum is always linear for suffi- 
ciently large q]. Termination can be viewed as a con- 
sequence of the restriction to positive sample measures 
/(a) > 

In the localized regime, the states contributing to the 
LDOS at a given position in the sample have a discrete 
energy spectrum, quantized by the typical localization 
volume £| 2 C . As a result, all non- unity LDOS moments 
diverge in the absence of level smearing. In a tunneling 
experiment, smearing can appear due to inelastic scatter- 
ing (temperature), open sample boundary conditions, or 
due to the finite energy resolution of the instrument. To 
characterize an Anderson insulating state over anLxL 
field of view with L 3> £| oc , the full LDOS distribution 
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should be cxamincd j 44 i 45 sensitive dependence of the dis- 
tribution shape to smearing can serve as a telltale sign of 
the localized regime. LDOS fluctuations in the Anderson 
insulator are reviewed in more detail in Sec. IIVBI 



B. Spin LDOS spectra 

By restricting the character of the tunneling species, it 
may be possible to measure individual LDOS components 
separately. For example, in the case of a spin-polarized 
(ferromagnetic) STM tip^ the spin-projected compo- 
nents z/f ^ can be separately resolved. The use of an un- 
polarized tip recovers the composite LDOS v = i/f + u^. 
We define the spin LDOS along the spin space direction 



i/(e,r) = i/f(e,r) - v[{e,v). 



(2.7) 



For a time-reversal invariant system (with or without 
spin-orbit scattering and/or disorder), one has v L {e, r) = 
0. In a system with broken time-reversal (e.g., magnetic 
impurities), but zero average spin polarization, the in- 
tegral of i/(e, r) over a sufficiently large region becomes 
arbitrarily small; we cannot use the normalized construc- 
tion in Eqs. (|2.1[) and (|2.2p to characterize spin LDOS 
multifractals. Instead, we employ the un- normalized in- 
verse spin participation ratio (ISPR) 

Pfe) = ]T (ni) q , £ = f d 2 r v\e, r). (2.8) 
In the scaling limit, 

Pfc)~«*(£)* S ~ 2 . (2-9) 



where the exponent x q is the scaling dimension for the 
corresponding moment operator in the disorder-averaged 
field theory description, and c q ^ for even q. 



III. WEAK DISORDER MULTIFRACTALITY 

A. Model. Short- and long-range correlated 
potential landscapes 

The Dirac surface states of a Z 2 topological insula- 
tor (TI) are guaranteed to appear in an odd number of 
flavorsi2i£ In this paper, we consider the simplest case of 
a single flavor, relevant to (e.g.) Bi 2 Sc3 and Bi 2 Tc3. The 
Hamiltonian is (in units such that H = 1) 



H = / d 2 rV> 1 



vf <5> [— idp + A M (r)] 
+ M(r) a 3 + V{r) 



1>, 



(3.1) 



where fi G {1,2}, and repeated indices are summed. The 
coordinates r = {x, y} chart the TI surface, while the 



topological bulk resides in the perpendicular z direction. 
In Eq. p.ip . Vp denotes the Fermi velocity, and the Dirac 
pscudospin Pauli matrices <r are related to the physical 
spin 1/2 operators S via {17^,(73} = 2{e fl „S u , S3}. The 
vector, scalar, and mass potentials {A, V, M} describe 
the effects of external electromagnetic fields and/or sur- 
face impurities. In the absence of time-reversal (7") sym- 
metry breaking, A = M = 0. (Sec Appendix [A] for an 
enumeration of discrete symmetry operations.) Thus, a 
mass gap is explicitly forbidden so long as T remains a 
good symmetry, a consequence of the protection afforded 
by the topologically nontrivial bulk. When T is broken 
by an external magnetic field B, the vector and mass 
potentials are 



A 



M 



-e4° rb) 



7|| 
2v F 



7-L 



B, 



(3.2) 



where 711 (xl) denotes the Zeeman coupling to the in- 
plane field By (out-of-plane field B z ), and the orbital 



effect is embedded in Aa'^ via e a pd a A^ 

Non-magnetic adatoms or charge traps are encoded in 
the scalar potential V(r). In-plane (out-of-plane) polar- 
ized magnetic impurities additionally induce point ex- 
change coupling to the vector A(r) [mass M(r)} fields^ 
The different types of disorder leading to V, A, and M 
are sketched in Fig.[TJ Assuming a random surface distri- 
bution of impurities and spatial rotational invariance on 
average, the disorder potentials can be taken as Gaussian 
white noise distributed variables, 



.(orb) 



B z 



V(r)V(v') = A v v F 6(r-r'), 
A a {r)Ap{v') = A A v\ 5 aP S(v - r'), (3.3) 



M(r)M(r') = Am vp 5(r — r'). 

The dimcnsionlcss variances Ay,A,M quantify the disor- 
der strength. In the first Born approximation, these are 
of the form 



At), 



(3.4) 



where rti mp is the impurity density, and u(q) denotes the 
Fourier transform of the single impurity potential. We 
note that a net in-plane magnetization of the surface im- 
purities A^ ^ can be removed by a gauge transforma- 
tion, while the average scalar potential V is absorbed into 
the chemical potential. We will assume that there is no 
net magnetization perpendicular to the surface, M = 0, 
or that we only probe LDOS fluctuations on energy scales 
much larger than the induced gap 2vpM . 

In 2D, the single impurity potential u(r) [Eq. (13.41) ] 
must decay faster than 1 jr 2 (or oscillate rapidly enough) 
so that the limit u(q — > 0) exists; otherwise, the white 
noise assumption in Eq. (|3.3[) is invalidated by long range 
impurity potential correlations 4£ This causes a problem 
for charged impurities, which can become poorly screened 
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for a small surface doping relative to the Dirac point. 
In graphene, the long-range correlated potential undu- 
lations induced by poorly-screened substrate impurities 
leads to a smearing of the Dirac point over an energy 
scale ksT rms cx vpy/n^, and to the breakup of the sam- 
ple into electron and hole puddles . 34 ' 35 The advent of 
electron-hole puddles has until recently prevented the ob- 
servation of various "intrinsic" phenomena associated to 
the Dirac carriers in graphene experiments such as veloc- 
ity renormalization 3 ^ and hydrodynamic transport near 
the Dirac point. In this respect, a large surface or bulk 
doping actually improves the situation for STM mea- 
surement of disorder-induced quantum interference, since 
these carriers screen the potential of surface charges. The 
disorder potential can be considered short-range corre- 
lated for scales larger than the screening length. 

If we consider only surface doping, with an insulating 
bulk, then the Thomas-Fermi wavelength due to a finite 
surface carrier density n is given by 

Atf = -J^-, (3.5) 
a V n 

where a = e 2 jevp is the effective "fine structure con- 
stant." The permittivity e = (1 + exi)/2, the average 
of the bulk TI below and vacuum above the surface. For 
Bi2Se3 with a surface density of n = 7x 10 12 cm -2 (corre- 
sponding to a doping level of 0.3 eV relative to the Dirac 
point), 30 v F = 5x 10 5 m/s (Ref.|30|), and permittivity 4 ! 
£ti = 113, one obtains Atf ~ 90 ran. This is very large, 
and indicates that the surface state carrier density is in- 
adequate to screen charged impurities. A smaller screen- 
ing length is possible for bulk doping ^ or by performing 
experiments on thin film samples exfoliated over a metal- 
lic gate. Alternatively, one can restrict the deposition of 
surface impurities to non-doping adatoms, e.g. iron in 
Bi 2 Se 3 ^ The disorder variance associated to Thomas- 
Fermi screened charged impurities is 

A V=7T 1 ^. (3.6) 
n 

Finally, we note that the appearance in isolation 
of any of the three disorder potentials in Eq. p.ip 
realizes three different symmetry classes of Anderson 
( derealization,— i 48 i 49 see Appendix [A] for a review. The 
T-invariant case with Aa.m = belongs to the spin-orbit 
class All, which is also the class of the Z 2 topological 
bulk [Fig. UK a)]. In the case of broken T, Ay, a/ = 
realizes the random vector potential model in class AIII 
[Fig.[T](b)], while Ay,A = gives the random mass model 
in class D [Fig. [T](c)]. All three classes exhibit delocal- 
ized states in 2D, although this occurs only at the Dirac 
point for class AIII.— In the T-invariant symplectic case, 
the unpaired single Dirac flavor avoids the usual spin- 
orbit metal-insulator transition,— remaining delocalized 
even for strong disorder due to a topological term^ The 
generic case of broken- T with all three disorder poten- 
tials non-zero realizes the unitary class A, and is believed 



to flow under renormalization to the plateau transition 
in the integer quantum Hall effecti 24 ' 29 (See Sec. IIV A 21 
for a review). 

Because in-plane (out-of-plane) Zeeman coupling ap- 
pears in the vector (mass) potential [Eq. (|3.2[) ], one is 
tempted to identify class AIII (class D) with the limit 
of an otherwise clean surface, dusted with charge neutral 
magnetic impurities randomly polarized in-plane (per- 
pendicular to the TI surface). However, a magnetic 
adatom is expected to also induce a local scalar potential 
deformation V(r). For example, it can dope the surface 
or bulk, as occurs for a manganese impurity in Bi2Sc3 
(Ref. |37T)] - As discussed in Appendix [3] the advent of 
any two flavors of disorder destroys the additional dis- 
crete symmetries enjoyed by the special class D and AIII 
Hamiltonians. The asymptotic long-distance LDOS scal- 
ing is then governed by the unitary class A, discussed 
above. Nevertheless, depending upon the relative mi- 
croscopic strength of the magnetic versus potential per- 
turbations induced by polarized magnetic impurities, the 
class AIII or D model may provide an adequate approxi- 
mation for broken- T LDOS fluctuations on intermediate 
scales. 



B. Results 

To compute the scaling of LDOS moments in a quan- 
tum theory with quenched disorder, one employs a path 
integral Z to express products of fermion Green's func- 
tions as functionals of the disorder configuration. Using 
a trick (replicas^ ' 20 ' 22 supersymmetry^ or Kcldysh 5 -^) 
to normalize Z = 1, the Green's functions are for- 
mally averaged over disorder configurations (typically 
with a Gaussian weight). The result is a translationally- 
invariant, but "interacting" field theory, where the disor- 
der strength A appears as a coupling constant £££ Per- 
turbative calculations are controlled via loop expansion 
for small A. 

To determine the scaling, one decomposes the q th 
LDOS moment into projections upon the renormaliza- 
tion group (RG) cigenoperators of the disorder-averaged 
theoryi22r— The multifractal spectrum r(q) is deter- 
mined by the most relevant (negative)^ scaling dimen- 
sion x q exhibited by an eigenoperator in this decomposi- 
tion, and is given b y 19 ' 43 

T{q) = 2(q-l)+x q -qx 1 . (3.7) 

1. Broken T: random vector potential disorder (Class AIII) 

The properties of the model in Eq. (|3.1[) with short- 
range correlated disorder [Eq. (|3.3[) ] were originally stud- 
ied in Ref. [24]- In this work, the exact multifractal spec- 
trum r(q) was calculated for the broken- T, random vec- 
tor potential (~ in-plane polarized magnetic impurity)^ 3 - 
class AIII model, to all orders in A a- Technically, this 
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result obtains because the disorder-averaged AIII model 
is conformally invariant at the Dirac point, and the ex- 
act LDOS moment spectra can be extracted using an 
Abelian bosonization treatment. The exact spectrum 2 ^ 
is quadratic in q, and takes the form of Eq. (|2.3[) . with 



0a = — ■ 

7T 



(3.8) 



Subsequent wor k 41 ' 42 on the random vector potential 
model elucidated the mechanisms of termination and 
freezing, transitions that occur in the spectral statis- 
tics for large moments q > q c (AA) or strong disorder 
A A > 2tt. 

For this broken-T class, we can also examine the spin 
LDOS fluctuations, utilizing the same nonperturbative 
bosonization treatment employed in Ref. [2J. The spin 
LDOS z/(e,r) taken along an axis I in spin space was 
defined by Eq. (|2 . T[) . Moment fluctuations are charac- 
terized by the inverse spin participation ratio (ISPR) in 
Eq. (|2.8[) . the scaling limit of which is controlled by the 
dimension x q that appears in Eq. (|2.9[) . The out-of-plane 

ISPR Vg(e) is associated to the "mass" fermion bilinear 
v 3 = ip^a^ip. For the random vector potential model, 
the most relevant contribution to V 3 (s) carries the same 
scaling dimension that gives the composite LDOS scaling 
in Eqs. (|2~3| and 



A 



A 2 
-1 ■ 



(3.9) 



The chiral components of the in-plane spin LDOS are the 
energy-resolved U(l) Dirac current operators 



= v 1 ± iv 2 = rp ji a±tp- 



(3.10) 



Moments of these are RG eigenoperators that receive no 
corrections. The scaling of the associated ISPR is gov- 
erned by the disorder-independent (tree level) exponent 

xf = q. (3.11) 



Eqs. (jcTS)) . and (|3~TT|) are exact results that hold 

to all orders in A^. 



2. Broken T: random mass disorder (Class D) 

In the rest of this section, we provide new results 
for the broken-T, random mass (~ out-of-planc polar- 
ized magnetic impurity)^ class D model, the T-invariant 
class All model, and the generic broken-T unitary class 
A model. For weak disorder, none of these are confor- 
mally invariant, and we resort to perturbation theory. In 
this section we summarize results; some technical aspects 
are sketched in Appendix [B] The results obtained below 
hold only for small Ay,M <C 1, wherein the disorder ap- 
pears as a weak marginal perturbation (at tree level) to 
the clean Dirac surface band structure. 



For the broken-T case of random mass disorder (with 
Ay = A^ = 0), one obtains quadratic multifractality at 
one loop, again governed by Eq. (|2.3j) . with 



'M 



Am 
2tt 



O (A 



M 



(3.12) 



Moments of the out-of-plane spin LDOS operator v 3 — 
ip^a^ip, as well as of the chiral in-plane [U(l) current] 
operators = ip^a±'ip constitute RG eigenoperators at 
one loop, with scaling dimensions 



A 



q * 2tt 
xf =q+0 (A M ) 



C(A M ), 



(3.13) 
(3.14) 



Note that the first correction in Eq. (|3.13[) is positive (and 
linear in q); this should be contrasted with the AIII case, 
Eq. (|3.9[) above. On general grounds, the anomalous scal- 
ing dimension associated to the q > 1 moment of the 
composite LDOS, or any projected component thereof, 
must appear with a negative sign. The reason is that 
this quantity is associated to a moment of a normalized 
probability distribution^^ through Eqs. (|2.1j) and (|2.2|) . 
For a quadratic r(q) spectrum, this leads in particular to 
9 > in Eq. (12.31) [consistent with a positive, real dis- 
order variance— c.f. Eqs. (|3~8)l . (|3~T2"T) . and (f3TT5l) ]. By 
contrast, the spin LDOS is defined as the difference be- 
tween two orthogonal projections [Eq. (|2.7[) ]; for this rea- 
son, the first disorder correction to the scaling dimension 
in Eq. (|3 . 13|) is not required to appear with a particular 
sign. 



3. Non-magnetic disorder ( Class All) 

In the T-invariant case of scalar potential disorder, it 
turns out that no local operator (without derivatives) ex- 
hibits multifractal scaling to first order in Ay- For Dirac 
fcrmions, this applies to both LDOS and energy-resolved 
current moments. Physically, the weak influence of non- 
magnetic disorder is due to interference mediated by the 
Dirac pseudospin (equivalent to physical spin 1/2 on the 
TI surface). The Dirac pseudospin is also responsible 
for the suppression of backscattering from a single non- 
magnetic impurity^ Technically, this result is derived 
by mapping the one-loop renormalization process of local 
operators to the action of a certain spin-1/2 Hamiltonian 
Hy , and identifying renormalization group eigenopera- 
tors with states that diagonalize Hy^ (see Appendix IB"]) . 
As a result, to lowest order one observes plane wave scal- 
ing in the LDOS IPR [Eq. (gl])]. The spin LDOS van- 
ishes exactly, due to T. 

The first non-trivial correction to the LDOS r(q) ap- 
pears at two loops. To this order, the spectrum is again 
quadratic as in Eq. (|2.3|) . A straight-forward but labori- 
ous calculation gives the coefficient in this equation, 



3A 2 

e v = ^y- + o(A 3 



(3.15) 
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FIG. 2: Quadratic multifractality for isolated disorder fla- 
vors. The Renyi dimension d q = 2 — 6 q [Eqs. (|2.3|) and 
(|3.16[) ] is plotted for the exact vector potential (AIII), one- 
loop mass (D), and two-loop scalar potential (symplectic All) 
results [Eqs. l[3~%)). If5^2"]l . and (pTT5)) ]. The disorder strength 
is A = 0.05 for each case. The broken time-reversal class D 
and AIII corrections appear at order Am,a, while the (much 
weaker) time-reversal invariant class All correction begins at 
order Av. 



Since Ay oc n lmp [Eqs. (|3.4[) or (|3.6p ], we find that the 
non-trivial multifractal scaling begins at second order 
in the impurity density. This is qualitatively weaker 
than any of the broken-T regimes, where the quadratic 
multifractality appears already at first order, Eqs. Q3.8P 
and (|3.12p . This distinction between T- invariant and 
T-broken surfaces is our primary result, and can be 
tested directly in STM experiments by varying the con- 
centration of deposited surface disorder. Although the 
T-invariant case is not conformally invariant (for a dis- 
cussion of renormalization effects, sec Sec. IIII C[ below), 
the multifractal r{q) and f(a) spectra depend only upon 
a single parameter, the variance Ay. Eq. (|3.15[) can 
be extended to higher loops, allowing ever more pre- 
cise tests against numerics or experimental data within 
the perturbatively accessible regime. The multifractal 
spectrum therefore provides a unique fingerprint for the 
time-reversal invariant Dirac surface state of the Z 2 topo- 
logical insulator, in the presence of weak but otherwise 
generic non-magnetic disorder. The opposite limit of 
strong disorder for the T-invariant case is discussed be- 
low in Sec. IIVA1I 



4- Broken T: generic disorder (Class A) 

When T is broken and any two disorder flavors appear, 
the system resides in the unitary class A. The third disor- 
der flavor is always generated under renormalization — see 
Sec. HITCl below. The results of Eqs. $3Ji§ and (|3~T2j) for 
the LDOS r(q) spectrum in the random vector and mass 
potential models suggest that the unitary case also ex- 
hibits multifractality to first order in the impurity density 

"-imp, Since A A ,M,V « "imp- 

With multiple flavors of the disorder, solving the op- 
erator mixing problem for the q th LDOS moment re- 



quires the diagonalization of an effective spin Hamilto- 
nian iT eff \ transcribed in Eq. (fBTOj) of Appendix [Bj In 
Figs. [3] and SI we present the results obtained by numer- 
ically diagonalizing this matrix for various combinations 
of {Ay, Am, A^i}. In these figures we plot the Renyi 
dimension^ d q , defined for q ^ 1 via 



d =1® 
a q — . 

q-l 



(3.16) 



Figs. [3] and 2] show that the generic broken-T case is 
multifractal at one loop, and easily distinguished from 
the two-loop T-invariant result, in the limit of weak dis- 
order. [Note that Fig. [4] indicates that the r(q) spectrum 
is not purely quadratic in this general case.] It should 
therefore be possible to precisely distinguish the broken- 
T and T-invariant spectra experimentally, by observing 
the dependence of the deviation 2 — d q on n- lmp . The 
single-disorder flavor results for comparable strengths are 
plotted in Fig. [5] for reference. 

For the multidisorder unitary model, the same RG 
cigenoperators dominate the scaling of composite v and 
out-of-plane spin is 3 LDOS moments. The dimension x q 



that determines the spin LDOS scaling via Eq. (|2.9[) also 
enters into the LDOS spectrum in Eq. (|3.7p . leading to 
Figs. [3] and |31 By contrast, moments of the chiral spin 
LDOS components v [Eq. (|3.10p ] remain eigenoperators 
that acquire no corrections at one loop, 



xf =q+ O(A a A ) 



(3.17) 



a,j3£ {A,M,V}. _ 

As reviewed in the subsequent Sec. IIV A 21 for M = 
0, the generic broken-T model is believed to flow to 
the critical state at the integer quantum Hall plateau 
transition ! 24 ! 29 This state exhibits strong multifractality 
that has been extensively studied in numericSi 18 i 19 ' 26 " — 
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FIG. 3: One-loop Renyi dimensions [Eq. (|3.16[) ] in the bro- 
ken T, multidisorder unitary class A case, for various disorder 
strength combinations. These results were obtained by nu- 
merically extracting the largest positive eigenvalue from the 
effective spin Hamiltonian H eff in Eq. (|B10p (restricting the 
search to operators invariant under spatial rotations and re- 
flections); see Appendix[B]for details. The two-loop result for 
the T-invariant case is shown for reference. 
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FIG. 4: The same as Fig. [3l but for different unitary class 
disorder strength combinations. In the data presented here, 
Av > Am, a- Regardless, the one- loop spectrum obtained for 
either Am, a > is stronger than the two loop T-invariant 
case with Am = Aa = 0. The latter is also shown for refer- 
ence. 



C. Renormalization effects 

As discussed at the beginning of the previous section, 
the disorder-averaged Dirac surface state theory used to 
compute LDOS multifractal spectra is an "interacting" 
field theory, wherein the disorder strengths Av,a,m ap- 
pear as coupling constants (c.f. Appendix IB1 . Because 
these parameters are dimensionless, at weak coupling the 
disorder constitutes a marginal perturbation of the clean 
Dirac band structure. The one-loop RG equations for 
these parameters are given b y 24 ' 54 



dAA 1 A A 



dA 



A I 



dl 

dA v 
dl 



1 

7T 
1 



{2A A - A M ) (A 



A I 



{2A A + A y ) (A 



A I 



-Ay) 

Ay) 



(3.18a) 
(3.18b) 
(3.18c) 



where I = logL denotes the log of the RG length scale 
(e.g., the system size). Energy e scales as 



dine 
dl 



(3.19) 



where the (scale-dependent) dynamic critical exponent is 



z = 1 + — (2A A + A M 

Z7T 



Ay) + O (A a Ap) , (3.20) 



a,pe{A,M,V}. 

In this section, we use Eqs. (|3.18p - (|3.20p to derive the 
dynamical scaling of the disorder parameters Av.a,m{^)\ 
here energy e is measured relative to the Dirac point, not 
the Fermi energy. (From the point-of-view of the disor- 
dered Dirac theory, a finite energy above the Dirac point 
constitutes a relevant perturbation.—) Using the results 
obtained in the previous section, we thereby determine 
the enhancement or suppression of LDOS multifractality 
approaching the Dirac point, due to renormalization. 



1. Broken T: random vector potential disorder (Class AIII) 

For the random vector potential model with Ay = 
A M = 0, Eq. (pTTB) implies 



dA A 
~df 



0. 



so that A a = A A (constant), where A^ is the "mi- 
croscopic" value derived from the randomly polarized 
in-plane magnetic impurity distribution.— This result in 
fact holds to all orders in A^j^i in this case, the theory 
describing LDOS fluctuations at the Dirac point is con- 
formally invariant. Multifractality is neither enhanced 
nor suppressed as one moves away from the Dirac point , 
defined as e = 0. However, for non-zero energies e ^ 0, 
in an infinite size sample all states are in fact localized^ 
The localization length diverges upon approaching the 
band center as £i c(e) ~ s~ x / z . with z = 1 + Aa/tt 
[Eq. (|3~2TJ|) ]. Eqs. (JO) and for r(q) hold on scales 
smaller than £| c(e)- 



2. Broken T: random mass disorder (Class D) 



For the random mass model with A 



v 



A, 



dA 



A I 



A 2 

M 



dl 



OA 



AI) 



so that the disorder is marginally irrelevant at 
weak coupling.— Integrating this equation and using 
Eqs. (|3.19j) and (|3.20[) . we can compute the scaling of Am 
with energy. At energy scale T, we define A M } = Am(T); 
then for the smaller energy scale e (relative to the Dirac 
point), we obtain the logarithmic suppression 



AmO < T) 



■ A (0) 

M 



o 



(AV) 



AI 



log ( - 



(3.21) 



This equation holds for |1 — e/T| <C 1. In the limit as 
e — > 0, the disorder strength vanishes as 



In 



O 



T 



A(O) 

AI 



A(O) 

AI 



In" 



A 



A I 



(3.22) 



For small A M \ Eq. (|3.22p applies only at very small en- 



ergies e < T exp( 
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3. Non-magnetic disorder (Class All) 

Now we consider the T-invariant model. The flow 
equation for Ay is 



A 2 



+ O (A= 



(3.23) 



In contrast to the random mass, the random scalar po- 
tential is a marginally relevant perturbation to the clean 
band structured Examining lower and lower energy 
scales approaching the Dirac point, one observes stronger 
effects of the disorder. In the asymptotic scaling limit 
wherein the impurity potential strength becomes "large" 
(Ay > 1), numerical^ results and analytical^ arguments 
imply that the disordered T-invariant Dirac theory renor- 
malizes into the "conventional" symplectic metal. The 
metal is distinguished from the Dirac theory by its non- 
zero (and non-critical) density of states at zero energy^ 
and by its r(q) spectrum.— ^ We discuss the strong cou- 
pling LDOS multifractality below in Sec. IIV A 11 If at 
energy T, Ay(T) = AS' <C 1, then for a somewhat 
smaller energy e we obtain the logarithmic enhancement 



Ay(e<T) 



.A (0) 
v 



O 



A(O) 



(3.24) 



Eq. (|3.24p implies that renormalization strengthens 
multifractality approaching the Dirac point e = 0, for 
the T-invariant case. We emphasize that this has nothing 
to do with weak (anti-)localization. The latter occurs in 
the diffusive metallic regime with /cp^mfp 3> 1, where l m f p 
denotes the elastic mean free path. The diffusive regime 
obtains at strong coupling^ near the Dirac point kp — > 
[Sec. IIV A 11 below]. The impurity strength renormaliza- 
tion in Eq. (j3.24p is a quantum effect deriving from the 
clean band structure, in the "near-ballistic" regime ^ 



4- Broken T : generic disorder (Class A) 

In the generic case of broken T, with multiple disor- 
der coupling strengths non-zero, the system flows toward 
strong coupling Av,u,A — > oo. As a result, multifractal- 
ity is enhanced approaching the Dirac point. The RG 
flow ultimately terminates at a strong coupling critical 
point, or in the Anderson insulator, discussed in the next 
section. 



IV. STRONG COUPLING REGIMES 

In this section, we review prior results on strong cou- 
pling regimes relevant to the disordered Dirac Z2 topo- 
logical insulator surface states, LDOS fluctuations and 
associated multifractal spectra. These are not new, but 



provide complimentary information to the new results 
derived in the previous section. 

In both generic cases of T-invariant, and T-breaking 
impurities, the disordered Dirac description used in 
Sec. IIIII fails on the largest length and lowest energy 
scales (approaching charge neutrality). For a sufficiently 
dilute concentration impurities, the results obtained in 
the previous section characterize the start of the scaling 
regime, over energy and length scales such that the dis- 
order strengths remain weak Av,a,m(L 7 e) <C 1. When 
these parameters become order one (due to renormal- 
ization down to lower energies and longer lengths), the 
system crosses over to one of the strong coupling regimes 
discussed below. 



A. Delocalized states at strong disorder 

1. T -invariant case: diffusive metal via strong disorder 

In a random scalar potential field, the Dirac point vac- 
illates in energy with spatial location; as a result, the den- 
sity of states near charge neutrality is enhanced by the 
disorder. Due to the suppression of pure backscattering 
for Dirac fermions^ the state density enhancement more 
than compensates for the increased scattering introduced 
by the additional impurities. As a result, scalar potential 
disorder actually increases the (zero temperature, Lan- 
dauer) conductance at charge neutrality beyond the clean 
ballistic result, e 2 /Trh (Refs. HIH). As in Sec. [ml here 
we assume short-range correlated disorder, due either to 
charge neutral impurities or efficient screening by bulk 
and/or surface carriers. We do not discuss the puddle 
regim o 34 ' 35 in the present paper. 

The effective disorder strength Ay is enhanced by 
renormalization, as indicated by the runaway flow im- 
plied by Eq. (j3.23p . The concomitant density of states 
and conductance growth suggests that the disordered 
Dirac theory ultimately crosses over to the ordinary dif- 
fusive symplectic metal, a result born out by numerics^ 
In the absence of time-reversal symmetry breaking, An- 
derson localization is prohibited on the surface of a topo- 
logical insulator 

The symplectic metal possesses a finite (non-critical) 
average density of states at charge neutrality, and a dis- 
tinct r(q) spectrum. For a large effective diffusion con- 
stant D (induced for a Dirac fermion subject to suffi- 
ciently strong disorder^ or for a weakly disordered sys- 
tem examined on large length scales), the lowest order 
result for the multifractal spectrum appears in Eqs. (|2 . 3[) 
and Eq. (|2.4j) . above. In the latter equation, j3 = 4 for 
the symplectic classi 23 ' 40 

For the T-invariant case, the strongest multifractal- 
ity is expected at intermediate coupling. Weak disorder 
Ay *C 1 induces weak multifractality in the Dirac lan- 
guage [r(q) in Eqs. (|2.3|) and (J3.15[) ]. while strong dis- 
order ultimately pushes the system into the symplectic 
metal, where a large diffusion constant D suppresses the 
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first correction in Eq. (|2.4 



for disorder-averaged LDOS moments evaluates to 



2. Broken T: IQHP transition 

For generic T-breaking disorder, i.e. all three Ay,M,A 
non-zero, the disordered Dirac theory is also unsta- 
ble under renormalization. When the average mass is 
zero M = (see below), the flow in Eq. f|3 . 1 8[) is be- 
lieved to terminate at the critical point of the integer 
quantum Hall plateau transition . 24 ! 29 This is the delo- 
calized state separating adjacent Hall plateaux; it ex- 
hibits strong multifractality that has been extensively 
studied in numerics . 18 ' 26 " — The spectrum is believed to 
be universal,— and is approximately^ parabolic as in 
Eqs. (1231) and ((231), with 9 ~ 0.26 (Refs. USEF 



B. Anderson insulator 

At zero chemical potential relative to the Dirac point, 
an average out-of-plane spin magnetization at the sur- 
face of a Z2 TI corresponds to the presence of a non- 
zero Dirac mass M for the surface carriers. This in- 
sulating state resides in a quantum Hall plateau [with 
o~xy = sgn(M) e 2 /2h]££&2^ In the presence of surface 
disorder, the plateau state will assume the character of 
a localized Anderson insulator. In this section we review 
LDOS fluctuations in the Anderson insulator. The dis- 
cussion is relevant not only to the magnetized surface of 
a 3D TI, but also to localized states populating the bulk 
gap of a disordered TI. Proposals exploiting localization 
to realize so-called "topological Anderson insulators" by 
adding impurities to clean hosts include those in Ref. l55l . 

To understand local density of states fluctuations in an 
Anderson insulator, it is useful to first study a toy prob- 
lem. Consider a tight-binding model on a d— dimensional 
lattice, subject to nearest-neighbor hopping t, and a ran- 
dom on-site potential Vi, distributed uniformly over the 
region —W/2 < Vi < W/2. We assume the absence of 
spatial correlations in the disorder potential. The inverse 
relative strength of the disorder is measured by the ratio 
t/W. We consider first the extreme limit of zero hopping, 
t/W — > 0. In that case, the LDOS is the on-site operator 



p(u) 



W/2 



W/2 



dV 

~w 



S[v-v{e,V)\ 



1 



TTV 2 W V V. 



(4.1) 



In this equation, the LDOS is constrained to the interval 
^min < v < is m3Xl where 



7r (W 2 + 4?7 2 ) 



1 

TTT] 



(4.2) 



Using Eq. (|4. 1 [) . one can compute the disorder-averaged 
moments of the LDOS, 



vi 



W-Ki-y^Tiq) 



n 



(4.3) 



The average LDOS is V = 1/W; all higher moments are 
proportional to rr x ^ q , and thus diverge in the limit of zero 
energy smearing r\ — > 0. This not surprising, because the 
energy spectrum in our trivial toy model is discrete, so 
that the LDOS operator becomes a delta function with 
ill-defined moments as r\ — > 0. The moments arc domi- 
nated by the power-law (Parcto) tail of the distribution, 
accumulating at the upper limit v = z^ max . By contrast, 
the typical LDOS, defined as v typ — cxp(logz^) is domi- 
nated by the infrared 



y ty P 



4?7e 2 
nW 2 ' 



This vanishes in the limit 77 — > 0. 

We see that observables exhibit broad statistics in the 



power- 



single site model, governed by the p(v) ~ v~^l 2 
law distribution in Eq. (|4.1[) . The moments are rendered 
finite only by the non-zero energy smearing parameter 
rj. This should be compared to the LDOS statistics in 
a system with extended states and weak multifractality, 
e.g. that characterized by the quadratic r(q) spectrum 
in Eq. (|2~3]l . with < < 1. It is knowr3 that the 
corresponding LDOS distribution has a Gaussian bulk, 
with a small amplitude log-normal tail responsible for 
the weak multifractality. For the metallic system, the 
result is independent of energy smearing, provided that 
the thermodynamic limit is taken before the smearing is 
set to zero. Returning to the toy insulator model, we 
observe that the global density vq of states (GDOS) is 
self-averaging in the same limit. The GDOS is defined 



r)/w 



{e-ViY + if 



1 N 



where 77 denotes an energy-smearing parameter. Physi- 
cally, smearing is determined by inelastic scattering, open 
sample boundary conditions, or due to the finite energy 
resolution of the probing instrument. 

At the "band center" e = 0, the distribution function 



where N denotes the number of sites. In the large ^V- 
limit, the cumulant expansion can be evaluated via the 
saddle-point. The cumulants of the GDOS then take the 
form 



L^gJ 



N 1 " 1 (7/? + ...). 
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where [■■•]* denotes the q th cumulant, and the omitted 
terms are smaller by positive powers of rj. Taking the 
infinite system size limit N — > oo before sending the en- 
ergy smearing to zero 77 — >- leads to the vanishing of all 
q > 1 GDOS cumulants. 

The calculations above can be extended to non-zero 
hopping via a locator expansion in small t/W, as per- 
formed by Anderson in his original 1958 paper^ This 
expansion can be formally summed to all orders in ID 
and on the Bethe lattice^ but an explicit solution for 
the LDOS statistics is difficult to obtain this way; see 
Ref. for an alternative approach. 

Altshuler and Prigodinss succeeded in deriving the dis- 
tribution generating disorder-averaged LDOS moments 
in a ID system, which is exponentially localized for ar- 
bitrarily weak disorder £ In the thermodynamic limit for 
a closed sample, they obtain the "inverse Gaussian" dis- 
tribution 




exp 



4tj {v - if 



(4.4) 



where v = u/V, and e is the typical energy level spacing 
in a localization volume; e _1 is also the clastic scattering 
lifetime^ 4 . In the limit of small smearing 77 <C e, this 
distribution has moments 



4 1 "«r (q - 



l-q 



(4.5) 



mcnts, one must employs* 



(loc) (?) = - 



rflnL 



In 



= / d d v V|V>i(r)| a «*(e- ei 

v J L d 



(4.6) 



Since the levels are discrete, 

ffiu> ?? )9- V(e,r) - V |^(r)| 29 <5 (e - *) . (4.7) 



Thus, 



(loc) (9) = - 



dlnL 

d 

d\nL 



In 



d d r 



lim(7r?7) 9 " 1 ^9(e,r) 

r)— >0 



lnjy d d r^(e,r)| 



(4.8) 



In this equation, we replace averagcs-of-the-logs with 
logs-of-the-average, a procedure that is legitimate here 
because spatial and disorder-averaging are expected to 
yield the same results on the insulating side. Noting that 
the LDOS moments are L-independent in the insulator 
for L ^> £| oc , we obtain the expected result^ 8 , for localized 
states 



r (loc) (9)=0, 
computed in a well-defined ry — ^ limit. 



(4.9) 



The exact result in Eq. (|4.5p for the ID Anderson in- 
sulator exhibits the same singular dependence on the 
energy smearing rj as the single site model moment in 
Eq. (|4.3| . In fact, the distributions in Eqs. (|4. 1[) and 
(|4.4[) are very similar: both feature a z^ -3 / 2 power law at 
intermediate v, while the exponential factor in Eq. (|4.4j) 
plays the role of the hard cutoffs ^min.max in Eqs. (|4.1[) 
and (|4.2[) . The close resemblance of the exact ID and 
single site model results can be attributed to the discrete 
spectrum of energy levels contributing to the LDOS in 
an Anderson insulator, with an energy level spacing de- 
termined by the localization volume ^ c in d spatial di- 
mensions. 

The take away is that the LDOS distribution in an 
Anderson insulator becomes very broad, with a power- 
law tail yielding divergent moments, in the limit of van- 
ishingly small energy smearing. In an STM experiment 
performed at ultra-low temperature, on a large, isolated 
Anderson localized surface, the collected LDOS statistics 
should be very sensitive to the smearing induced by the 
energy resolution of the measurement itself. 

The locally discrete energy spectrum of the LDOS in 
the Anderson insulator invalidates the use of Eq. (|2.2[) 
as a tool to compute the multifractal r(q) spectrum. As 
advocated above, the shape of the LDOS distribution 
function and its sensitivity to smearing can best reveal 
the insulating phase. If one insists upon computing mo- 
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Appendix A: Discrete symmetries, random matrix 
classification, and disorder 

The 10 symmetry classes of disordered Hamiltoni- 
ans (Hermitian random matrices) can be efficiently dis- 
tinguished by the presence or absence of time-reversal 
T, particle-hole V, and chiral/"sublattice" symmetry 
Ci 7 ' 48 ' 49 For the two-component Dirac Hamiltonian in 
Eq. (|3.1[) . the definitions of these symmetries are es- 
sentially unique. In terms of the two-component Dirac 
spinor if), these appear as 

T : ij} — > —iaixj), i — > —i (Ala) 

V: V-^iiy] 1 ", ( Alb ) 

C: ip -> 63 [i/' t ] , (Ale) 



12 



In the second quantized language, T and C are antiuni- 
tary transformations; the unitary V can be taken as the 
product of these. 

The imposition of any one of the discrete symmetries 
upon the Hamiltonian in Eq. (|3.1j) in every disorder real- 
ization restricts its form, and selects a particular ran- 
dom matrix symmetry classi 7 ' 24 ' 48 ' 49 (1) T-invariancc: 
A = M = 0, only potential disorder Ay > is al- 
lowed. Since T 2 = — 1, this is the symplectic (spin- 
orbit) class All, which is also the symmetry class of 
the (presumed T-invariant) topological Z2 bulk. (2) V- 
invariance: V = A = 0, only random mass disorder 
Am > is allowed. Because V 2 = +1, this is the bro- 
ken time-reversal class D. (3) C-invariance: V = M = 0, 
only random vector potential disorder A^ > is allowed. 
This is the broken time-reversal class AIII. (Technically, 
it is the "topological" /WZW class AIIIi in the language 
of Refs.0®) 

Class All is generically realized whenever time-reversal 
is unbroken. Magnetic impurities randomly polarized 
parallel (perpendicular) to the TI surface manifest as 
point exchange sources in the vector (mass) potentials 
of Eq. (|3.1[) ; we are thus tempted to identify symme- 
try classes D and AIII with these two limits. How- 
ever, a magnetic impurity will typically induce a local 
potential fluctuation V(r) as well. As a consequence, 
the generic case of broken-time reversal symmetry corre- 
sponds to the absence of T, V, and C, which gives the uni- 
tary class A.— i 48 i 49 In fact, for a vanishing average mass 
M = 0, the surface of a topological insulator with generic 
time-reversal breaking disorder is expected to flow under 
renormalization to the critical point of the integer quan- 
tum Hall plateau transition ! 24 ' 29 On a different note, the 
class AIIIi and class D versions of H in Eq. (|3.1|) can be 
realized on the surface of a bulk T-invariant 3D topolog- 
ical superconductor, where time-reversal is respectively 
preserved or broken at the surface.— 



Appendix B: Perturbation theory 



1. Chiral Decomposition and one-loop 
renormalization 



We write a 2+0-D fermion path integral to represent 
correlation functions in the disordered Dirac Hamilto- 
nian transcribed in Eq. (|3.1|) . The fermion operators are 
replaced with the Grassmann fields {tp,^} — > {V'ijV'i}; 
here i £ {1, . . . , n} denotes a replica index, and we are to 
send n — ► at the end of the calculation.^ 2 ^ We employ 
a "chiral decomposition" of the two-component spinors, 



ipi = [B4 Li] 



(Bl) 



Then the action of the replicated theory is 

' e (R^i + L l R i ) + RiLicj) + LiRi 



S 



tf r 



Ri {-id + A)R, + L l (-id + A) L t 



(B2) 



where we have introduced complex coordinates {z, z} = 
x ± iy, {8,8} = (d x =F id y ), and disorder potentials 
{A, A} = A x =F iA y , {(j), <f>} = V ± M, The energy e is a 
fixed parameter. In Eq. (|B2[) . repeated replica indices are 
summed. Assuming the Gaussian white-noise variances 
for the disorder potentials enumerated in Eq. (|3.3|) . the 
replicated theory can be averaged over disorder configu- 
rations. The post-ensemble averaged action is 



S = Sq+Sa+S\ + S2, 
where Sq is the clean Dirac action, and 

S A = - 2A a J (frRiRiliLi 
Ay + A M 



Si 



So 



(Ay - Am) 



err (RiLiRjLj 



d v Ri Li Lj Rj . 



(B3) 



(B4a) 



LiRiLjRj) , 



(B4b) 
(B4c) 



Different replicas become coupled through the 
disorder 

The disorder-averaged composite LDOS v(e, r) corre- 
sponds to the fermion bilinear expectation 



77= (ipip) = (RL + LR). 



(B5) 



The spin LDOS v~ u (e, r) was defined by Eq. (|2~7)) . For the 
out-of-plane and in-plane (chiral) components, one has 

=$& 3 i>) = (RL - LR), (B6a) 
1A = 2{(AR>, (LL)}, (B6b) 

The overlines appearing in the left-hand sides of 
Eqs. (|B5I) and (|B6|) denote disorder-averaging, whereas 
the angle brackets on the right-hand sides represent in- 
tegration in the fermion path integral, using the action 
S in Eq. (lB3l) . 

A generic local operator corresponding to the q th mo- 
ment of some fermion bilinear can be viewed as sum of 
"strings," where each string consists of 2q total right (R) 
and left (L) mover labels, arranged in some order. For 
example, the disorder-averaged q th moment of the LDOS 
is represented by the the composite operator expectation 
value 



I* (r) = (lllRiLiir) + LiRi(r)] 



(B7) 



In this equation, a product is taken over operators carry- 
ing indices in the first q < n replicas. The q th LDOS mo- 
ment is computed by placing one copy of the LDOS op- 
erator into each of q different replicas; before averaging, 
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this gives v q in a fixed realization of the disorder. (Plac- 
ing instead the q copies into the same replica would give 
the disorder-averaged first moment of a 25-point Green's 
function.) The operator in Eq. (|F37|) is an even weight 
sum of 2 q "strings", all of length 2q. I.e., 



i/9 (r) ={RL; RL; . . . ; RL} 

+ {LR; RL; RL;...; RL} 
+ {RL; LR; RL; . . . ; RL} 
+ ... 

+ {LR;LR;...;LR}. 



(B8) 



The semicolons separate fermion bilinears in different 
replicas. Each bilinear has two entries, corresponding to 
the chiral identity of the barred and unbarred operators. 

The set of all length 2q strings forms a complete ba- 
sis for q th moment local operators (without derivatives). 
These basis strings mix under renormalization due to 
the disorder^ In general, the composite operator (= 
weighted string sum) corresponding to the q th moment 
of a bilinear as in Eq. (|B7[) does not constitute an eigen- 
operator of the renormalization group. The main task is 
to (1) identify RG cigcnopcrators for each disorder type 
and compute the spectrum of scaling dimensions, and (2) 
compute the projection of the LDOS and (for broken T) 
spin LDOS moment operators onto this eigenbasis, and 
determine the most relevant contributions. 



a. Effective Hamiltonian for 1-loop renormalization 

It is useful to view each string BjS ci configuration of 
2q spin 1/2 moments. We associate {R, R} — > 1/2 (spin 
up) and {L, L} — > —1/2 (spin down). Operators invari- 
ant under spatial rotations have equal numbers of up and 
down spins, and therefore reside in the zero total magne- 
tization sector with J z = 0. We picture each string as a 
basis state for a length q chain, with two spins per site. 
Sites are labeled by the replica index i 6 {1, . . . , q}. The 
two spins at each site are distinguished by labels "A" 
and "B," corresponding to barred and unbarred opera- 
tors, respectively. 

Renormalization occurs via the action of the disorder 
vertices appearing in Eq. (|B4[) , employing the clean Dirac 
propagator in a standard loop expansion. Operator mix- 
ing at one-loop is encoded in the effective "Hamiltonian" 



#(eff) 



In A 
""2V 



^2 ( S Ai ~ S Bi) ( S Aj ~ S 



Bj) 



(A M + A y ) £ (5+ Sb S + S+ Bl S 



i,j=l 
1 



(A M -A v ) £ (5+5^+5+5 



7 Bj 



(B9) 



In this equation, S% Bi denotes a spin-1/2 operator act- 
ing on the barred (A) or unbarred (B) spin in replica i. 
The prefactor obtains from evaluating the loop integrals 
using a hard momentum cutoff A. The first, second, and 
third lines in the heavy brackets arise through the action 
of the disorder vertices in Sa, Si, and S2 , respectively. 
The A a renormalization is diagonal in the 'f/i (R/L) ba- 
sis. By contrast, Si, and £2 perform single exchanges of 
right and left labels. Si (S2) mediates interflavor AhB 
(intraflavor A <H> A, B <H- B) exchanges. Summing the 
angular momenta, 



In A 

~~2n 



2A A [J z a - J B ) 

+ 2(A M + A V )(J%J% + JIJ%) 
+ (A M - A y ) 



(J1Y 



(J z bY 



(BIO) 



where 3 a ,b = J2t &A,Bi- 

In the general case of broken T discussed in Sec. MI B~4l 
all three disorder parameters arc present. The most rele- 
vant eigenvalue of Eq. (|B10[) determines the scaling of the 
g th LDOS moment^ The first few multifractal moments 
for various disorder configurations were obtained through 
numerical diagonalization; results appear in Figs. [3] and 

m , - 

Even moments of the out-of-planc spin LDOS u 3 

[Eq. (|B6al) ] arc invariant under spatial rotations and 
parity^ In the multidisorder unitary case, even mo- 
ments of the composite v and out-of-plane spin v 3 LDOS 
are dominated by the same RG eigenoperator. The di- 
mension x 3 that determines the spin LDOS scaling via 
Eq. (|2.9[) is the same that enters into the LDOS spec- 
trum in Eq. (|3.7[) . Figs. |3]and|31 Moments of the chiral 
spin LDOS components in Eq. (|B6b[) correspond to the 
highest weight states \j = q;m = ±q); here, j(j + 1) de- 
notes the eigenvalue of (3a + 3b) 2 , with < j < q and 
—j < m < j. These highest weight states are annihilated 
by in Eq. ()B"T0|) . leading to Eq. (|3~T7) . 

Below we discuss the special cases of isolated disorder 
flavors. 



b. Broken T: random vector potential disorder ( Class AIII) 



For Am = Ay = 0, Eq. (|B10|) reduces to 



H 



(eff) 



lllA A I T z T z \ 2 
A A (J A ~ Jb) 

7T 



In A . 
= A A (m A 

IT 



m B ) 



(Bll) 



On the second line, we have evaluated H^f f ' for the 
product state \ja, Jb] m A, tub)- Since maxfj^.g) = q/2 
and |j71a/b| < Ja/Bi the maximum eigenvalue attains 
for the states \q/2, q/2; q/2, -q/2) -> {RL; RL;...; RL} 
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and \q/2,q/2;-q/2,q/2) -> {LR; LR; . . . ; LR}. These 
have J z = 0, and thus correspond to operators invari- 
ant under spatial rotations; the symmetric combination 
is also parity- invariant^ Via standard renormalization 
group machinery ) 60 ' 61 one obtains the most relevant scal- 
ing dimension for a q— fold product operator, 



4 A) =q - q 2 



A, 



(B12) 



Using Eq. (|B12j) in Eq. ([577]) gives the result for the 
quadratic r(g) spectrum in Eqs. (|2.3|) and (13. 



c. Broken T: random mass disorder (Class D) 



The Hamiltonian in Eq. (|B10[) with Am 
comes 



(eff) 



Ay In A 
2tt 



J 2 



zn2 



(J*) 



be- 



(B15) 



where J = + J^. Eq. (|B15|) has the same form as 
Eq. (|B13[) . with Ay — > —Am- This is consistent with a 
mapping between the random mass and vector potential 
models identified in Ref. H^- In the case of the scalar 
potential, the maximum eigenvalue is associated to the 
highly degenerate singlet sector \j = 0,m = 0), leading 
to the scaling dimension 



For the random mass case, Eq. (|B10j) becomes 



H 



(eff) 

M 



Am In A 
: 2^ 

A M In A 
' "hi 



i) 



(B13) 



where J = 3a + Jb- On the second line, we have eval- 
uated the "Hamiltonian" on a total angular momentum 
cigenstate \jm). For 2q spins, we have max(j) = q. The 
maximum eigenvalue is associated to the non-degenerate 
j = q, m = state, which is invariant under spatial 
rotations. The scaling dimension is 



r (M) 



2tt ' 



(B14) 



The corresponding eigenoperator \j = q,m = 0) is an 
equal weight symmetric sum of all permutations of q 
"i?" and q "i" labels, and has non-zero overlap with 
the naive LDOS moment (a g-fold triplet product) in 
Eq. (jB7)) . Using Eq. (|BT4)) in Eq. (JXTJ), one obtains the 
result for the quadratic r(q) LDOS spectrum in Eqs. (|2.3[) 
and (|3.12l) . By contrast, the out-of-plane spin LDOS 
(mass operator) in Eq. (|B6ap is a singlet; the disorder- 
averaged q th moment thus corresponds to the eigenoper- 
ator \j = 0, m = 0) [leading to Eq. (f3~T3)) ]. 



d. Non-magnetic disorder (Class All) 

We rotate the "B" composite spin by ir around the 
£-axis, 

TX , TX — JX JV . TV — _ TV JZ , JZ JZ 

J B > J B — J B , J B r J B — J B , J B t J B = J B . 



2tt 



(B16) 



Using this result in Eq. (|3.7|) gives r(q) = 2(q — 1). 
We conclude that no moment operator (without deriva- 
tives) acquires multifractal scaling at one loop for the 
T-invariant model. 



2. Two loop renormalization, T-invariant class All 

In the T-invariant class All model, the first correc- 
tion to the LDOS r(g) spectrum appears at second order 
in Ay. We have carried out a two- loop calculation and 
found that the naive LDOS moment in Eq. (|B7|) remains 
an eigenoperator. To this order, we find the scaling di- 
mension 



A 2 
^y 



2tt 8tt 2 



Ml - 1) 



O 



(A 3 



v)- 
(B17) 



Combining Eqs. fl3J]) and (|Bl7|) . we recover the 
quadratic multifractality for r(q) quoted in the text, 
Eqs. (|2.3[) and (|3.15|) . We have used dimensional rcgular- 
ization to obtain the result in Eq. (|B17|) . Although the 
Clifford algebra becomes formally infinite 6 - 2 - upon dimen- 
sional continuation to d = 2 — e, this causes no problems 
for the renormalization of the q th LDOS moment because 
the latter is already an eigenoperator. We omit details 
of the (lengthy) two-loop calculation in this paper. 
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